#!/usr/bin/env python
""" 
This function generates the TranSEQ input from
(1) targetGene.isoforms (Name and length for each enumerated isoforms)
(2) targetGene.mappings (Mapping to those isoforms)
"""
from one_AceView_gene_one_isoforms_file import *
from get_unique_AceviewGene_mappings import *
from PerM_output_format import *
from mapping2TranSEQ import MappingsOfIsoforms2AGene
import commands
import platform
import os

def get_read_length( mapping_file ):
    ''' Given a mapping file from PerM, return the read length
    of the first mapping.'''
    DEFAULT_READ_LENGTH = 100
    if( os.path.exists( mapping_file ) ):
        ifile = open( mapping_file , mode = "r" )
        line = ifile.readline()
        if ( line.strip() == "" ):
            return( DEFAULT_READ_LENGTH )
        mapping = PerMMappingLineInfo( line )
        read_length = len( mapping.read )
        if read_length > 0:
            return( read_length )
    return( DEFAULT_READ_LENGTH )


def get_all_mapping4aGene_from_all_mapping_files( gene_name ):

    gene_mapping_filename = gene_name + ".mapping"
    if( os.path.exists( gene_mapping_filename ) ):
        os.remove( gene_mapping_filename )
    if platform.system() == 'Windows':
        # for test purpose only
        mapping_files = [ 'm0.mapping', 'm1.mapping' ]
    else:
        cmd_result = commands.getstatusoutput( "ls *.mapping" )
        mapping_files = cmd_result[1].split()
    for file in mapping_files:
        get_mapping4aGene( file, gene_name )

def get_Isoforms( gene_name ):
    isofroms = []
    ifile = open( gene_name + ".isoforms", mode = "r" )
    for line in ifile:
        info = line.split()
        isofroms.append( ( info[0], int( info[1] ) ) )
    return( dict( isofroms ) )

def get_TranSEQ_input_4_a_gene( gene_name ):
    '''
    This function integrate three steps to generate the input for TranSEQ
    # (0) get_all_mapping4aGene_from_all_mapping_files()
    #     create the isoform_name and length file from fasta 
    (1) Get the (isoform name, length) dictionary from a pre-generated file
    (2) Get the reads mapped to the target gene from all mapping file 
    (3) output the input for TranSEQ (.dat)
    '''
    # get_all_mapping4aGene_from_all_mapping_files( gene_name ) (0)
    # create the isoform_name and length file from fasta
    isoform_lengthes_dict = get_Isoforms( gene_name ) # (1)    
    mappings = MappingsOfIsoforms2AGene( gene_name )
    mappring_filename = gene_name + ".mapping"
    ifile = open( mappring_filename , mode = "r" )
    for line in ifile:
        if ( line.strip() == "" ):
            break
        mapping = PerMMappingLineInfo( line ) # (2)        
        try:
            geneName = extract_AceView_gene_name( mapping.isoform_name )
            if geneName == mapping.isoform_name : #if not AceView format
                # Assume all mappings are from the same gene
                geneName = mappings.geneName
            mappings.add_mapping( geneName, mapping.read_name, mapping.isoform_name )
        except AttributeError:
            print line + " has no attribute "
            break
    # update the read length because the output matrix depend on the read length
    mappings.read_length = get_read_length( mappring_filename )
    mappings.get_TranSEQ_input( isoform_lengthes_dict ) # (3)
    ifile.close()

def usage():
    sys.stderr.write( "Working dir is " + os.getcwd() + '\n' )
    print __doc__
    print "The only arg is the gene name that will be analyzed for TranSEQ."

def main():
    try:
        opts, args = getopt.getopt( sys.argv[1:], 'h', ["help"] )
    except getopt.error, msg:
        print msg
        sys.exit( 2 );
    # process options and argument 
    for o in opts:
        if o in ( "-h", "--help" ):
            print __doc__
            sys.exit();
    if( len( args ) == 1 ):
        get_TranSEQ_input_4_a_gene( args[0] )
    else:
        usage()

if __name__ == '__main__':
    main()
